LONG WAVE RUNUP ON RANDOM BEACHES 
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Abstract. The estimation of the maximum wave rmiup height is a problem of practical 
importance. Most of the analytical and numerical studies are limited to a constant slope 
plain shore and to the classical Nonlinear Shallow Water (NSW) equations. However, 
in nature the shore is characterized by some roughness. In order to take into account 

(^ I the effects of the bottom rugosity various ad-hoc friction terms are usually used. In this 

Cn ■ paper we study the effect of the roughness of the bottom on the maximum runup height. A 

pi i' stochastic model is proposed to describe the bottom irregularity and its effect is quantified 

D i using Monte-Carlo simulations. For the discretization of the NSW equations we employ 

^^ I modern finite volume schemes. Moreover, the results of the random bottom model are 

f^ ■ compared with the more conventional approaches. 
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1. Introduction 



The estimation of the long wave runup on a sloping beach is a practical problem which 
O '• attracts nowadays a lot of attention due in part to the intensive human activity in coastal 

^^ I areas. The main demand comes from the coastal and civil engineering but also from coastal 

communities which are exposed to the tsunami wave hazard [SB06]. Consequently, a lot 
of effort is devoted to the development of fast and accurate estimation methods of the 
wave runup and horizontal excursion over a sloping beach [TS96, KS06, DP08, MSIO]. In 
general this problem is solved in simplified geometries (e.g. constant slope beach) and in 
the framework of Linear or Nonlinear Shallow Water (LSW, NSW) equations. However, 
^ . more general situations may require the application of other models and different numerical 

" ■ techniques (see e.g. [LWL02, MBFS07, DPDIO, DKMll] and the references therein). 

In practice, the available data are always subject to some uncertainties. For example, 
the bathymetry is known only in a discrete number of scattered points, while in reality 
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the shores are characterized by some rugosity. The missing information can be modeled 
by the inclusion of random effects. These circumstances have lead several authors to 
consider water wave propagation in random media [GJP93, dBCDE+OS, NaclO]. In the 
present study we model the natural beach roughness by small random perturbations of 
the smooth average bottom profile. The long wave dynamics are described by the classical 
NSW equations. We note that the dispersive effects could also be included (see [DKMll]), 
however they do not modify qualitatively the results that follows bellow. The main effect 
of the dispersion is a small reduction of the maximum runup height due to the wave energy 
flux to shorter wavelengths. 

2. Mathematical model and results 

Consider an incompressible perfect fluid layer bounded below by the solid bottom d{x) 
and above by the free surface ri{x,t). In the present study we are interested in the long 
wave regime which is described by the NSW equations: 

Ht + {Hu)., = 0, (2.1) 

{Hu)t+{Hu^ + lH^)^ = gHd^-gHSf, (2.2) 

where H{x,t) = d{x) + ri{x,t) is the total water depth and u{x,t) is the depth-averaged 
fluid velocity. The channel bottom d{x) is assumed to be a sloping beach described by the 
depth function d{x) = do — tan 6 ■ {x + £), where S is the constant bottom slope and i is the 
half-length of the physical domain. Parameters do, i, S are chosen so that a dry sloping 
area is below the still water level (see Table 1). The term Sf is included to model some 
friction effects and it will be taken zero unless otherwise noted. We consider the Boundary 
Value Problem (BVP) posed on the one-dimensional interval X = [—i,i], where on the 
right boundary x = £ we impose the so-called wall boundary condition u{i,t) = (in our 
simulations the wave front does not achieve this point), while on the left end x = — ^ we 
generate an incoming wave of height ri{—L,t) = —ao sm{ijJot)'H(To — t), where Tiit) is the 
Heaviside step-function and Tq = 27i/uo is the wave period. In other words, we generate 
a shorewared traveling, one-period monochromatic leading depression wave. The values of 
the various physical and numerical parameters used in this study are given in Table 1. 
The interval X is divided into cells Ci = [Xj„i ,a;j_|_i] of length Axi = Xi_^_i — x^_i and 

Xi = |(a;j_i + x^_^_l) denotes the midpoint of Cj, i = 1,...,N. Without any loss of 

generality we assume the partition of cells T = {Ci}fLi is uniform. In order to model the 
bottom roughness we construct a random perturbation in the following way. Let us fix 
an integer number r > 1 which will be referred to as the regularity parameter. Then, on 
each cell {Cjr}Y=i C T, with m = [yj, we generate a normally distributed pseudorandom 
variable C,j ~ A/'(0,o"^), where the parameter a characterizes the perturbation magnitude 
since \C,j\ < 1.96 • cr with probability 95%. Constructed in this way the random vector 
$, = {^j}Y=i is interpolated on the whole grid using cubic splines, for example, to obtain 
a particular realization of micro irregularities. The discrete bathymetry function becomes 
di = do — tan 5 ■ (xj + £) + C,i on each cell Cj. Several realizations of the random bottom for 
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Table 1 . Various parameters used in this study. Note that with the present 
choice of parameters do and g we solve the governing equations (2.1), (2.2) 
in the dimensionless form. 



Parameter 


Value 


Domain half-length, L 


17.0 


Bottom slope, tan 5 


0.06 


Gravity acceleration, g 


1.0 


Water depth at the left end, do 


1.0 


Incoming wave amplitude, oq 


0.15 


Incoming monochromatic wave frequency, Uo 


0.2 


Number of control volumes, N 


1000 


Number of Monte-Carlo runs, M 


1000 



Random sloping beaches 




Figure 1. A sample realization of random bottoms for a = 5 x 10~^ and 
various values of the regularity parameter r = 1, 2, 4, 5 and 8 starting 
correspondingly from the top. 



various values of r are shown on Fig. 1. If r = 1 we obtain a white noise while increasing 
this parameter is equivalent to the application of a spectral filtering operation (see Fig. 2). 
The hyperbolic system of NSW equations is discretized using the finite volume method, 
cf. [DKMll]. Specifically we use the characteristic flux approach, [GKCOl], combined with 
the UN02 space reconstruction procedure, [H087]. The well-balancing of the scheme is 
achieved by applying the well known hydrostatic reconstruction method [AB05]. The run- 
up algorithm description can be found in [DPDIO, DKMll]. For the time discretization 



D. DUTYKH, C. LABART, AND D. MITSOTAKIS 



Spectral density 




50 100 150 200 250 300 350 400 450 500 

k 



Figure 2. Spectral density of the bottom perturbation for several values of 
the regularity parameter r = 1, 2, 8 (starting from the top). 



we use the 3rd order Bogacki-Shampine Runge-Kutta scheme with adaptive time step 
selection. 

Once the parameters a and r have been chosen, we can generate a particular realization 
of the rough sloping beach and solve the BVP to determine the maximum wave runup. 
The shoreline motion of one particular realization with a = 10~^ and r = 1 is represented 
on Fig. 3. For comparison, the shoreline behavior in the idealized smooth bottom case is 
also represented on Fig. 3 with the blue dashed line. One can see that the main effect of 
the bottom rugosity is the reduction of the maximum wave run-up height -Rmax- In this 
particular simulation the wave run-up has been reduced by a factor of 2 approximately. 
Sometimes this effect is referred to as the apparent diffusion, cf. [NaclO]. Intuitively we 
can understand this outcome since a wave dissipates more energy due to the interaction 
with these micro irregularities. 

One of the main questions we address here is to quantify the run-up reduction when the 
bottom roughness varies. Our approach consists of performing direct numerical simulations 
of this process over random bottoms instead of adding some ad-hoc terms to model this 
roughness. We will return to this point below. In probabilistic terms we would like to 
estimate the expectation E,{R^ax) over all possible realizations of the random bottom noise. 

Since a random bottom perturbation is constructed in the discrete space, the dimension 
of the random parameters vector ^ G M™" scales with the number of control volumes N 
in our spatial discretization of the interval X. More precisely m = [A^/rJ , where r > 1 is 
the noise regularity parameter introduced above. The discrete space in our simulation is 
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Figure 3. The shoreline motion in the case of a smooth shore (the blue 
dashed line) and a particular realization of the random bottom with a = 
10~^, r = 1 (black solid line). 



of dimension A^ which is typically of the order 10^ (see Table 1). This value is imposed 
by the accuracy requirements of our direct simulations and this rather high dimension is 
a limiting factor for the choice of the expectation E(i?niax) numerical method estimation. 
Popular nowadays the polynomial chaos expansion method does not apply if the number of 
random parameters is typically greater than two. The Quasi-Monte- Carlo approach fails for 
dimensions higher than 200 because of substantial difficulties to generate a low discrepancy 
sequence of random vectors of such a large dimension. Consequently, we are limited to the 
standard Monte-Carlo method which is not sensitive to the stochastic problem dimension. 
However, we can apply a variance reduction method described below. 

In order to estimate E(i?max)5 we simulate M random bottom realizations, and for each 
case j we compute numerically the maximum run-up -RmL. We approximate E(i?max) by 
the mean Sm '■= jj 'Yl,j=i -RmL. According to the central limit theorem we know that 



E(-Rmax) belongs to the interval \S 



confidence. 



where a\i 




+ 1.96i 



with a 95% level of 



1 \^M /pi 

•- A/-1 2^^=11-"-" 

To reduce the size of the confidence interval, 



is an unbiased converging estimator of 
the variance of -Rmax- To reduce the size of the confidence interval, we can either increase 
M (which requires more computational time) or try to find a random variable with mean 
lE(-Rmax) and variance smaller than Var(_Rinax)- We opt for the second possibility - the 
so-called variance reduction technique. Since -Rmax can be seen as a function R{$,), where $, 
follows a centered Gaussian law ^/{Om, cr'^^m), we can use the adaptive importance sampling 
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Maximum runup value 
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Figure 4. Maximum runup value as a function of the perturbation char- 
acteristic magnitude a for the irregular case r = 1 (blue line) and the reg- 
ularized noise r = 6 (black line). For comparison, the red dash-dotted line 
represents the maximum runup value for the smooth bottom case. 



technique proposed in [LLll]. This method uses the fact that: V0 G M"*, E(i?(^)) 



E(i?(^ + e)e 



-e-i-i^. 



Then, one can construct an algorithm which finds the parameter 



vector 0* minimizing the variance of H{0,^) := R{^ + 6)e ^'^' 



\e\^ 



Then, the average 



value E(i?inax) is approximated by Sm '■= jj X]j=i -^(^i-i; ^j); where {Oj}Y=i is a sequence 
converging to 0* . We refer to [LLll, Section 2.2] for theoretical results on the central limit 
theorem in this adaptive case where the random variables are not independent anymore. 
This algorithm allows us to reduce the variance by a factor of two approximately. In our 
computations the confidence interval length has never exceeded 0.5% of the corresponding 
maximum run-up value with parameter M specified in Table 1. The probability density 
function of the -Rmax distribution for r = 1 and two values of a (10~^ and 10~^) are depicted 
on Figure 6. 

The Monte-Carlo simulation results are presented on Figs. 4 and 5. The dependence of 
the maximum run-up -Rmax value on the roughness magnitude a for two fixed values of the 
noise regularity r = 1 and 6 is shown on Fig. 4. On the other hand, the dependence of 
-Rmax on the regularity parameter r for several fixed values of a is represented on Fig. 5. 
We can see that the bottom roughness reduces significantly the wave run-up height while 
the noise regularization has an antagonistic effect. 
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Maximum runup value 




Figure 5. The maximum runup value as a function of the regularity param- 
eter r for several values of the perturbation magnitude a (increasing from 
the top). 



Since stochastic Monte-Carlo simulations of the bottom rugosity are computationally 
expensive, various friction ad-hoc terms are used to model these effects. The following 
examples can be routinely found in the literature: 

Chezy law:: Sf = c/^, where c/ is the Chezy friction coefficient 

-g^, where A is the resistance value determined ac- 



Darcy-Weisbach law:: Sj 

cording to the Colerbrook- White relation: 4= 

Manning-Strickler law:: Sf 

cient. 



^2 u\u\ 



-2.031og(^,,3,^^ 
where c/ is the Manning roughness coeffi- 



The friction coefficient c/ measures the bottom roughness as the parameter a in our random 
bottom roughness construction. Consequently, we can ask the same question: how does the 
maximum run-up value depends on the friction coefficient c/ if this term is incorporated 
into the model? We perform a series of deterministic numerical simulations for various 
values of c/ and the maximum wave run-up -Rmax is being measured. The numerical results 
are presented on Fig. 7. We can see that Chezy and Darcy-Weisbach laws provide a strong 
friction which reduces considerably the maximum run-up height. However, the Manning- 
Strickler law shows qualitatively a very similar behaviour to the results predicted by our 
stochastic model in the non-regularized case r = 1. 
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Figure 6. Probability density distribution for cr = 10 ^ (left image, 
E(i?max) = 0.78088, 5% and 95% quantiles are equal to 0.77734 and 0.78416 
respectively) and a = 10"^ (right image, E(i?max) = 0.42281, 5% and 95% 
quantiles are equal to 0.37857 and 0.46261 respectively). Note the difference 
in the vertical scales and in the horizontal extent of two distribution func- 
tions. 



3. Conclusions 

In the present study we considered the long wave run-up problem over rough bottoms. 
Specifically, we proposed a stochastic model to mimic the natural bottom roughness. Using 
the Monte-Carlo variance reduction technique, we quantified the maximum wave run-up 
behavior for various practically important values of the noise magnitude and regularity a, r. 
The maximum run-up is monotonically decreasing as the bottom roughness parameter a 
increases. However, this apparent dissipative effect might be drastically reduced when the 
noise regularity r is increased. Namely, in our simulations we observed the difference of the 
factor about two between the maximum run-up on the irregular (r = 1) and regularized 
(r = 6) perturbations. These results indicate that the regularity parameter has to be taken 
into account in some way while designing coastal protecting structures. Since the recent 
field survey by Fritz et al. [FKM"'"07] it has been known, for example, that coastal forests 
do not provide effective damping to tsunamis. 

Moreover, our stochastic computations were compared to several simulations using clas- 
sical friction terms routinely used to model the bottom rugosity. A very good qualitative 
agreement (for r = 1) was obtained with the Manning-Strickler law, while the Chezy and 
Darcy-Weisbach laws provide too strong momentum damping. 
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Figure 7. Comparison of the runup reduction effect for various ad- hoc fric- 
tion terms and the random bottom perturbation model in the non-regularized 
case r = 1. The horizontal axis represents the friction coefficient C/ for de- 
terministic computations and a for the random roughness model (blue solid 
line) . 
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